Each Genome

Concept

Load human genome split of samples, labeling with hs1..hs4. Output matrix.
Convert to common mouse genes. Re-load into object. Load mouse genome split of samples, labeling with ms1..ms4. Subset overlapping genes. Merge all 8. Normalize and scale. Show clusters labeled by cell; split by sample. Show nCount_RNA also. Determine significantly-different list mouse vs human

hs1.data=Read10X(data.dir="./hg19/S1/outs/filtered_feature_bc_matrix/")
hs2.data=Read10X(data.dir="./hg19/S2/outs/filtered_feature_bc_matrix/")
hs3.data=Read10X(data.dir="./hg19/S3/outs/filtered_feature_bc_matrix/")
hs4.data=Read10X(data.dir="./hg19/S4/outs/filtered_feature_bc_matrix/")
ms1.data=Read10X(data.dir="./mm10/S1/outs/filtered_feature_bc_matrix/")
ms2.data=Read10X(data.dir="./mm10/S2/outs/filtered_feature_bc_matrix/")
ms3.data=Read10X(data.dir="./mm10/S3/outs/filtered_feature_bc_matrix/")
ms4.data=Read10X(data.dir="./mm10/S4/outs/filtered_feature_bc_matrix/")

#some human gene symbols have underscores (but these are not in geneTrans).
#Substitute a dot so as not to raise an error.
rownames(hs1.data)=gsub("_",".",rownames(hs1.data))
rownames(hs2.data)=gsub("_",".",rownames(hs2.data))
rownames(hs3.data)=gsub("_",".",rownames(hs3.data))
rownames(hs4.data)=gsub("_",".",rownames(hs4.data))

#rename cells
colnames(x=hs1.data) <- paste('hs1',colnames(x=hs1.data),sep="_")
colnames(x=hs2.data) <- paste('hs2',colnames(x=hs2.data),sep="_")
colnames(x=hs3.data) <- paste('hs3',colnames(x=hs3.data),sep="_")
colnames(x=hs4.data) <- paste('hs4',colnames(x=hs4.data),sep="_")
colnames(x=ms1.data) <- paste('ms1',colnames(x=ms1.data),sep="_")
colnames(x=ms2.data) <- paste('ms2',colnames(x=ms2.data),sep="_")
colnames(x=ms3.data) <- paste('ms3',colnames(x=ms3.data),sep="_")
colnames(x=ms4.data) <- paste('ms4',colnames(x=ms4.data),sep="_")

#create objects
hs1=CreateSeuratObject(counts=hs1.data,project="MG",min.cells=5)
hs2=CreateSeuratObject(counts=hs2.data,project="MG",min.cells=5)
hs3=CreateSeuratObject(counts=hs3.data,project="MG",min.cells=5)
hs4=CreateSeuratObject(counts=hs4.data,project="MG",min.cells=5)
ms1=CreateSeuratObject(counts=ms1.data,project="MG",min.cells=5)
ms2=CreateSeuratObject(counts=ms2.data,project="MG",min.cells=5)
ms3=CreateSeuratObject(counts=ms3.data,project="MG",min.cells=5)
ms4=CreateSeuratObject(counts=ms4.data,project="MG",min.cells=5)

#Use https://satijalab.org/seurat/essential_commands.html 
hg=merge(x=hs1,y=c(hs2,hs3,hs4),project="Hg")
mg=merge(x=ms1,y=c(ms2,ms3,ms4),project="Mg")

#at this point we don't need all the original sample objects--can re-create if needed
rm(hs1,hs2,hs3,hs4,hs1.data,hs2.data,hs3.data,hs4.data)
rm(ms1,ms2,ms3,ms4,ms1.data,ms2.data,ms3.data,ms4.data)

#summary of each object
print(hg)

## An object of class Seurat 
## 13283 features across 16708 samples within 1 assay 
## Active assay: RNA (13283 features)

print(mg)

## An object of class Seurat 
## 17304 features across 29055 samples within 1 assay 
## Active assay: RNA (17304 features)

Convert human genes to mouse

Output raw data to data frame. Use conversion table (geneTrans) to translate human symbols to mouse symbols.

#load gene translation table
geneTrans=read.table("geneTrans.txt",sep=",",header=T,stringsAsFactors = F,row.names = 1)

mito.m=grep("^mt",rownames(mg),value=T)
mito.h=grep("^MT-",rownames(hg),value=T)
#turns out, these two vectors (n=13 each) are ordered and match, build an extended translation table
mitoTrans=data.frame(row.names = paste("mm10",mito.m,sep="_"),
                     Human.Symbol = mito.h,
                     Homologene_ID = rep(NA,13),
                     None = rep("yes",13),
                     Mouse.Symbol = mito.m,
                     hg19 = paste("hg19",mito.h,sep="_"),
                     mm10 = paste("mm10",mito.m,sep="_")
                     )

#extended translation table
xTrans = rbind(geneTrans,mitoTrans)

#extract human raw counts into table
hg.raw=GetAssayData(hg,slot="counts")

#subset rows in geneTrans (not necessary but simpler)
hg.raw=hg.raw[row.names(hg.raw) %in% xTrans$Human.Symbol,]  #cut from 13283 to 11364 rows

#translate human symbols to mouse
hg.trans=merge(x=hg.raw,y=xTrans[,c(1,4)],by.x=0,by.y="Human.Symbol",all.x=T)
rownames(hg.trans)=hg.trans$Mouse.Symbol
hg.trans=hg.trans[,!(names(hg.trans) %in% c("Row.names","Mouse.Symbol"))]

#convert matrix back into Seurat object
hm=CreateSeuratObject(hg.trans,project="MG",min.cells=5)

#summarize converted matrix
print(hm)

## An object of class Seurat 
## 11377 features across 16708 samples within 1 assay 
## Active assay: RNA (11377 features)

Subset mouse

Not all mouse gene symbols are available in the homologene table so filter table rows for this.

mg.raw=GetAssayData(mg,slot="counts")
mg.raw=mg.raw[row.names(mg.raw) %in% xTrans$Mouse.Symbol,]
mg=CreateSeuratObject(mg.raw,project="MG",min.cells=5)
print(mg)

## An object of class Seurat 
## 14545 features across 29055 samples within 1 assay 
## Active assay: RNA (14545 features)

Merge

#merge with mouse
hm=merge(x=hm,y=mg,project="HM")

#clean up objects no longer needed
rm(hg,mg,hg.raw,hg.trans,mito.m,mito.h,mg.raw)

#summarize merged object
print(hm)

## An object of class Seurat 
## 14781 features across 45763 samples within 1 assay 
## Active assay: RNA (14781 features)

Scale and normalize

Use standard workflow to calculate percent.mito (now all mouse symbols) and visualize by sample.

#filter and normalize merged object
mito.features=grep(pattern="^mt",x=rownames(x=hm),value=T)

percent.mito=Matrix::colSums(x=GetAssayData(object=hm,slot="counts")[mito.features,]) / Matrix::colSums(x=GetAssayData(object=hm,slot='counts'))
hm[['percent.mito']] = percent.mito

Diagnostic plots.

VlnPlot(object=hm,features=c("nFeature_RNA","nCount_RNA","percent.mito"),ncol=3)

FeatureScatter(object=hm,feature1 = "nCount_RNA",feature2 = "percent.mito")

FeatureScatter(object=hm,feature1 = "nCount_RNA",feature2 = "nFeature_RNA")

Scale

Print summary data before and after subsetting. Then normalize, find variable genes, and scale.

print(hm)

## An object of class Seurat 
## 14781 features across 45763 samples within 1 assay 
## Active assay: RNA (14781 features)

hm=subset(hm,nFeature_RNA > 200 & nFeature_RNA < 2500 & percent.mito < 0.05) #try > 100 & < 2500 & < .5
print(hm)

## An object of class Seurat 
## 14781 features across 17412 samples within 1 assay 
## Active assay: RNA (14781 features)

hm=NormalizeData(hm,normalization.method = "LogNormalize",scale.factor=1e4)
hm=FindVariableFeatures(hm,selection.method = 'mean.var.plot',mean.cutoff = c(0.0125,3),dispersion.cutoff = c(0.5,Inf)) #finds 4047 features
length(VariableFeatures(hm))

## [1] 3719

hm=ScaleData(hm,features=rownames(hm),vars.to.regress = c("nCount_RNA","percent.mito"))  #this takes a while

## Regressing out nCount_RNA, percent.mito

## Scaling data matrix

Dimension reduction

Start with PCA and determine how many dimensions are informative.

hm=RunPCA(hm,features=VariableFeatures(hm),verbose=T)

## PC_ 1 
## Positive:  B2m, Cd74, C1qb, Rpl13a, C1qc, C1qa, Spp1, Rpl29, Apoc1, H2-Ea-ps 
##     Cybb, Cx3cr1, C3, Fcgr4, Ifi203, Sat1, Ier2, H2-Eb1, Ms4a6b, Dusp1 
##     Ccl3, Angptl7, Ms4a7, Fos, Glipr1, Folr2, Junb, Rgs1, Ch25h, Lilrb4a 
## Negative:  Meg3, Syt1, Snap25, Epha5, Meis2, Nrxn3, Grin2b, Ahi1, Celf4, Atp1b1 
##     Arpp21, Scg5, Ndrg4, Rtn1, Snhg11, Il1rapl1, Pcp4, Plppr4, Pcsk2, Grm5 
##     Gpm6a, Camk2b, Synpr, Gad1, Gria3, Ly6h, Myt1l, Xist, Caly, Cadps 
## PC_ 2 
## Positive:  Cd74, C1qb, C1qc, C1qa, Spp1, Apoc1, Cybb, H2-Ea-ps, Cx3cr1, C3 
##     Rpl13a, Fcgr4, Sat1, Rpl29, Ms4a6b, H2-Eb1, Ccl3, Angptl7, Ms4a7, Glipr1 
##     Folr2, Lilrb4a, Adam28, Rgs1, Mis18bp1, P2ry12, Tgoln2, Clec7a, Ptafr, Ctss 
## Negative:  Flt1, Cldn5, Igfbp7, Itm2a, Id1, Bsg, Pglyrp1, Ptprb, Abcb1a, Ablim1 
##     Crip1, Adgrl4, Egfl7, Spock2, Pecam1, Epas1, Klf2, Ramp2, Wfdc1, Slc2a1 
##     Cxcl12, Slc7a1, Cyyr1, Hmcn1, Ahnak, Pltp, Foxq1, Prom1, Sox17, Sema3c 
## PC_ 3 
## Positive:  Syt1, Snap25, Meg3, Celf4, Tmsb10, Gad1, Ndrg4, Camk2b, Grin2b, Atp2b1 
##     Snhg11, B2m, Cd74, Bcl11a, Rpl13a, Myt1l, Slc12a5, Pcp4, Nrxn3, Rpl29 
##     Spp1, Ano3, Plppr4, Synpr, Gad2, Syt6, Sipa1l1, Atp1a3, C1qtnf4, Gria3 
## Negative:  Atp1a2, Slc1a2, Gja1, Car2, Clu, Cryab, Ptn, Qk, Cnp, Id4 
##     Cldn11, Cxcl14, Mobp, Sox9, Ermn, Acsl3, Rorb, Slc7a10, Mog, Aqp4 
##     Hes5, Tspan2, Slc6a11, Ppp1r14a, Slc4a4, Sept4, Opalin, Mag, Scrg1, Ugt8a 
## PC_ 4 
## Positive:  Gja1, Slc1a2, Atp1a2, Id4, Sox9, Clu, Cxcl14, Slc4a4, Rorb, Slc7a10 
##     Chchd10, Gpm6a, Tmem47, Aqp4, Fjx1, Hes5, Slc6a11, Vegfa, Acsl3, Lhx2 
##     Rcn2, Slc6a1, Sox2, Lix1, Cdh2, Rfx4, Itih3, Ezr, Id3, Oaf 
## Negative:  Mobp, Cldn11, Ermn, Mog, Tspan2, Tubb4a, Cryab, Sept4, Stmn4, Mal 
##     Ugt8a, Mag, Ppp1r14a, Kctd13, Opalin, Cnp, Tmeff2, Nkx6-2, Kcna1, Dixdc1 
##     Ptprd, Efnb3, Stmn1, Epb41l3, Edil3, Ppp1r16b, Mbp, Slain1, Anln, Aspa 
## PC_ 5 
## Positive:  Spp1, Cd74, C3, H2-Ea-ps, Cybb, Ifi203, Fcgr4, Apoc1, Tpm1, Adam28 
##     Folr2, H2-Eb1, Lilrb4a, Glipr1, Rpl29, Tgoln2, Clec7a, Mis18bp1, Ch25h, Ptn 
##     Olig1, Naip1, Vsig4, Ms4a7, Rnase1, Ms4a6b, Gimap4, Igf1, Epn2, Rpl13a 
## Negative:  Hexb, Ctss, Fyb, Ctsd, Cd52, P2ry6, Cx3cr1, C1qa, C1qc, P2ry12 
##     Rgs2, Ly6e, C1qb, Cd72, Eef1g, Ccrl2, Ier5, Tnfrsf13b, Junb, Il1b 
##     Wdr89, Atf3, Col27a1, Btg2, Slamf9, Icam1, Nfkbiz, Xist, Serpine1, Tanc2

print(hm[['pca']],dims=1:5,nfeatures=5,projected=F)

## PC_ 1 
## Positive:  B2m, Cd74, C1qb, Rpl13a, C1qc 
## Negative:  Meg3, Syt1, Snap25, Epha5, Meis2 
## PC_ 2 
## Positive:  Cd74, C1qb, C1qc, C1qa, Spp1 
## Negative:  Flt1, Cldn5, Igfbp7, Itm2a, Id1 
## PC_ 3 
## Positive:  Syt1, Snap25, Meg3, Celf4, Tmsb10 
## Negative:  Atp1a2, Slc1a2, Gja1, Car2, Clu 
## PC_ 4 
## Positive:  Gja1, Slc1a2, Atp1a2, Id4, Sox9 
## Negative:  Mobp, Cldn11, Ermn, Mog, Tspan2 
## PC_ 5 
## Positive:  Spp1, Cd74, C3, H2-Ea-ps, Cybb 
## Negative:  Hexb, Ctss, Fyb, Ctsd, Cd52

VizDimLoadings(hm,dims=1:2)

DimPlot(hm)

hm=ProjectDim(hm)

## PC_ 1 
## Positive:  B2m, Tyrobp, Cd74, C1qb, Tmsb4x, Rpl13a, Aif1, Ftl1, C1qc, C1qa 
##     Gpx1, Spp1, Cyba, Rpl29, Lst1, Apoc1, H2-Ea-ps, Cybb, Srgn, Laptm5 
## Negative:  Gria2, Meg3, Syt1, Pcsk1n, Snap25, Epha5, Meis2, Nrxn3, Grin2b, Ahi1 
##     Dclk1, Celf4, Atp1b1, Arpp21, Scg5, Stmn3, Dlgap1, Ndrg4, Rtn1, Snhg11 
## PC_ 2 
## Positive:  Tyrobp, Cd74, C1qb, C1qc, C1qa, Aif1, Spp1, Apoc1, Lst1, Cybb 
##     H2-Ea-ps, Cx3cr1, C3, Laptm5, Rpl13a, Fcgr4, Gpr34, Alox5ap, Sat1, Fcer1g 
## Negative:  Flt1, Cldn5, Igfbp7, Itm2a, Id1, Bsg, Pglyrp1, Ptprb, Abcb1a, Adgrf5 
##     Ablim1, Crip1, Adgrl4, Egfl7, Sox18, Jcad, Slc9a3r2, Esam, Spock2, Pecam1 
## PC_ 3 
## Positive:  Syt1, Mef2c, Snap25, Meg3, Celf4, Tmsb10, Gad1, Ndrg4, Camk2b, Grin2b 
##     Atp2b1, Snhg11, B2m, Cd74, Bcl11a, Pcp4l1, Rpl13a, Myt1l, Srgn, Slc12a5 
## Negative:  Plpp3, Atp1a2, Aldoc, Slc1a2, Gja1, Tsc22d4, Car2, Prdx6, Csrp1, Ntsr2 
##     Clu, Cryab, Gpr37l1, Plp1, Mt1, Slc1a3, Ptn, Qk, Mfge8, F3 
## PC_ 4 
## Positive:  Gja1, Slc1a2, Mt3, Mt2, Atp1a2, Id4, Aldoc, Plpp3, Slc1a3, Ntsr2 
##     Sox9, Clu, Cxcl14, Sparcl1, Prdx6, S1pr1, Atp1b2, Ptprz1, Mfge8, Slc4a4 
## Negative:  Plp1, Mobp, Cldn11, Ermn, Mog, Tspan2, Tubb4a, Cryab, Sept4, Stmn4 
##     Mal, Ugt8a, Mag, Ppp1r14a, Kctd13, Opalin, Tmem88b, Cnp, Tmeff2, Ttll7 
## PC_ 5 
## Positive:  Spp1, Cd74, C3, H2-Ea-ps, Neat1, Cybb, Ifi203, A2m, Fcgr4, Apoc1 
##     Tpm1, Arl5a, H2-DMb1, Adam28, Folr2, H2-Eb1, Lilrb4a, Olr1, Fcgrt, Glipr1 
## Negative:  Hexb, Ctss, Fyb, Ctsd, Cd52, Rnase4, Vsir, Rpl17, Lgmn, Rps17 
##     Tmem119, Fcgr3, Ctsz, Mafb, Selplg, Cst3, Ptpn18, Serinc3, Unc93b1, Rack1

DimHeatmap(hm,dims=1,cells=500,balanced=T)

DimHeatmap(hm,dims=1:6,cells=500,balanced=T)

hm=JackStraw(hm,num.replicate=100)
hm=ScoreJackStraw(hm,dims=1:20)
JackStrawPlot(hm,dims=1:20)

## Warning: Removed 52060 rows containing missing values (geom_point).

ElbowPlot(hm)

Start clustering

#start clustering
hm=FindNeighbors(hm,dims=1:13) #adjust dims based on plots

## Computing nearest neighbor graph

## Computing SNN

hm=FindClusters(hm,resolution=0.2)

## Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
## 
## Number of nodes: 17412
## Number of edges: 634880
## 
## Running Louvain algorithm...
## Maximum modularity in 10 random starts: 0.9751
## Number of communities: 15
## Elapsed time: 5 seconds

table(Idents(hm))

## 
##    0    1    2    3    4    5    6    7    8    9   10   11   12   13   14 
## 3760 1886 1667 1553 1551 1319  996  985  978  977  536  477  281  235  211

Find tSNE projection and add to object

hm=RunTSNE(hm,dims=1:13)
DimPlot(hm,reduction='tsne')

DimPlot(hm,reduction='tsne',split.by='orig.ident')

FeaturePlot(hm,features='Spp1')

FeaturePlot(hm,features='Hexb')

FeaturePlot(hm,features='nCount_RNA')

UMAP

hm=RunUMAP(hm,dims = 1:13)
DimPlot(hm,reduction='umap')

FeaturePlot(hm,features='Spp1')

FeaturePlot(hm,features='Hexb')

#FeaturePlot(hm,features='nCount_RNA',split.by='orig.ident')
FeaturePlot(hm,features='nCount_RNA')

Label clusters

Use top 2 genes from prior clustering to do this, following

hm.markers=FindAllMarkers(hm,only.pos=T,min.pct = .25,logfc.threshold = .25)

## Calculating cluster 0

## Calculating cluster 1

## Calculating cluster 2

## Calculating cluster 3

## Calculating cluster 4

## Calculating cluster 5

## Calculating cluster 6

## Calculating cluster 7

## Calculating cluster 8

## Calculating cluster 9

## Calculating cluster 10

## Calculating cluster 11

## Calculating cluster 12

## Calculating cluster 13

## Calculating cluster 14

hm.markers %>% group_by(cluster) %>% top_n(10,avg_logFC)

## # A tibble: 150 x 7
## # Groups:   cluster [15]
##    p_val avg_logFC pct.1 pct.2 p_val_adj cluster gene  
##    <dbl>     <dbl> <dbl> <dbl>     <dbl> <fct>   <chr> 
##  1     0      2.76 0.901 0.164         0 0       Slc1a2
##  2     0      2.53 0.945 0.218         0 0       Mt2   
##  3     0      2.45 0.882 0.093         0 0       Gja1  
##  4     0      2.41 0.612 0.271         0 0       Mt1   
##  5     0      2.39 0.583 0.068         0 0       Aldoc 
##  6     0      2.30 0.693 0.062         0 0       Id4   
##  7     0      2.29 0.947 0.265         0 0       Mt3   
##  8     0      2.27 0.573 0.06          0 0       Plpp3 
##  9     0      2.22 0.541 0.051         0 0       Clu   
## 10     0      2.20 0.643 0.155         0 0       Slc1a3
## # … with 140 more rows

Save top 5 list

hm.markers %>% group_by(cluster) %>% top_n(5,avg_logFC) %>% write.csv("Top 5 per cluster.csv",row.names = F)
hm.markers %>% group_by(cluster) %>% top_n(5,avg_logFC) %>% as.data.frame %>% paged_table

Identify human and mouse microglial clusters

Use Spp1 and Hexb to identify cluster number for human and mouse microglia, respectively.

humClus=as.character(subset(hm.markers,gene=="Spp1")$cluster)
print(humClus)

## [1] "1"

musClus=as.character(subset(hm.markers,gene=="Hexb")$cluster)
print(musClus)

## [1] "3" "5"

Find genes differentially expressed between human and mouse microglia

diffGenes=FindMarkers(hm,ident.1=humClus,ident.2=musClus,min.pct=0.1)

Display summary of differentiall expressed genes

table(sig=diffGenes$p_val_adj <= 0.05, twofold=abs(diffGenes$avg_logFC) >= 1)

##        twofold
## sig     FALSE TRUE
##   FALSE    64    0
##   TRUE   1033  161

diffGenes %>% tibble::rownames_to_column() %>% filter(p_val_adj <= 0.05) %>% filter(abs(avg_logFC) >= 1) %>% paged_table

Save differential expressed list as object (for later use in R) and csv (as text to preserve gene names)

saveRDS(diffGenes,"each_genome_diffGenes.rds")
write.table(as.data.frame(diffGenes),"each_genome_diffGenes.csv.txt",sep=",",quote=T)

Cluster Averages

hm.mean=AverageExpression(hm)

## Finished averaging RNA for cluster 0

## Finished averaging RNA for cluster 1

## Finished averaging RNA for cluster 2

## Finished averaging RNA for cluster 3

## Finished averaging RNA for cluster 4

## Finished averaging RNA for cluster 5

## Finished averaging RNA for cluster 6

## Finished averaging RNA for cluster 7

## Finished averaging RNA for cluster 8

## Finished averaging RNA for cluster 9

## Finished averaging RNA for cluster 10

## Finished averaging RNA for cluster 11

## Finished averaging RNA for cluster 12

## Finished averaging RNA for cluster 13

## Finished averaging RNA for cluster 14

head(hm.mean$RNA)

##                0          1          2           3          4           5
## A1bg  0.00000000 0.56872724 0.00000000 0.000000000 0.00000000 0.014836597
## A2m   0.15314837 5.18379849 0.04210029 0.005251912 0.09368136 0.297655861
## Aaas  0.05733664 0.11813280 0.08439040 0.055573357 0.14291236 0.036201081
## Aacs  0.12539327 0.03850084 0.01330111 0.048624545 0.00000000 0.007303951
## Aagab 0.05812589 0.07979784 0.00000000 0.282635428 0.01157240 0.002539016
## Aak1  0.41231548 1.95849839 0.16050586 0.304748043 0.30731988 0.112739540
##                6          7          8          9         10         11
## A1bg  0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000
## A2m   0.11449824 0.01356341 0.04508776 0.04600083 0.02872112 0.08915362
## Aaas  0.02594271 0.09654678 0.01881447 0.08820378 0.09924200 0.06282685
## Aacs  0.01148760 0.03355285 0.04207683 0.02500704 0.04811068 0.06270119
## Aagab 0.03326553 0.12360453 0.11109761 0.04593593 0.12846729 0.02146108
## Aak1  0.07711254 0.51745374 1.39913976 0.44800418 0.35380163 0.55236947
##                12         13         14
## A1bg  0.025750498 0.00000000 0.00000000
## A2m   0.473156957 1.83015150 0.00000000
## Aaas  0.073334448 0.05353564 0.02591217
## Aacs  0.008036854 0.00000000 0.00000000
## Aagab 0.030291339 0.05914882 0.02830128
## Aak1  0.257230746 0.07358431 0.34979784

write.csv(hm.mean$RNA,"hm cluster averages.csv.txt")

Show species by original split

#stash idents 
hm[["old.ident"]]=Idents(hm)
#get vector of cell idents
all.cells=Cells(hm)
#split by species
hg.cells=grep("^h",all.cells,value=T)
mm.cells=grep("^m",all.cells,value=T)
#apply new idents
Idents(hm,cells=hg.cells)="Human"
Idents(hm,cells=mm.cells)="Mouse"
table(hm$orig.ident,Idents(hm))

##      
##       Mouse Human
##   hs1     0  2419
##   hs2     0  1463
##   hs3     0  2633
##   hs4     0  2912
##   ms1  1787     0
##   ms2   776     0
##   ms3  3459     0
##   ms4  1963     0

Plot tsne labeled by species no legend

DimPlot(hm,label=T,repel=T,label.size=8,reduction = "tsne")+NoLegend()

Plot nCounts in tsne labeled by species

FeaturePlot(hm,features='nCount_RNA',reduction = "tsne")

FeaturePlot(hm,features='nCount_RNA',reduction = "tsne",split.by="ident")

Session Info

sessionInfo()

## R version 3.5.3 (2019-03-11)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.2 LTS
## 
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/atlas/libblas.so.3.10.3
## LAPACK: /usr/lib/x86_64-linux-gnu/atlas/liblapack.so.3.10.3
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] rmarkdown_1.12    dplyr_0.8.0.1     Seurat_3.0.0.9000
## 
## loaded via a namespace (and not attached):
##  [1] nlme_3.1-137        tsne_0.1-3          bitops_1.0-6       
##  [4] RColorBrewer_1.1-2  httr_1.4.0          tools_3.5.3        
##  [7] utf8_1.1.4          R6_2.4.0            irlba_2.3.3        
## [10] KernSmooth_2.23-15  lazyeval_0.2.2      colorspace_1.4-1   
## [13] npsurv_0.4-0        withr_2.1.2         tidyselect_0.2.5   
## [16] compiler_3.5.3      cli_1.1.0           plotly_4.9.0       
## [19] labeling_0.3        caTools_1.17.1.2    scales_1.0.0       
## [22] lmtest_0.9-36       ggridges_0.5.1      pbapply_1.4-0      
## [25] stringr_1.4.0       digest_0.6.18       R.utils_2.8.0      
## [28] pkgconfig_2.0.2     htmltools_0.3.6     bibtex_0.4.2       
## [31] htmlwidgets_1.3     rlang_0.3.4         zoo_1.8-5          
## [34] jsonlite_1.6        ica_1.0-2           gtools_3.8.1       
## [37] R.oo_1.22.0         magrittr_1.5        Matrix_1.2-17      
## [40] Rcpp_1.0.1          munsell_0.5.0       fansi_0.4.0        
## [43] ape_5.3             reticulate_1.12     R.methodsS3_1.7.1  
## [46] stringi_1.4.3       yaml_2.2.0          gbRd_0.4-11        
## [49] MASS_7.3-51.1       gplots_3.0.1.1      Rtsne_0.15         
## [52] plyr_1.8.4          grid_3.5.3          parallel_3.5.3     
## [55] gdata_2.18.0        listenv_0.7.0       ggrepel_0.8.0      
## [58] crayon_1.3.4        lattice_0.20-38     cowplot_0.9.4      
## [61] splines_3.5.3       SDMTools_1.1-221    knitr_1.22         
## [64] pillar_1.3.1        igraph_1.2.4        future.apply_1.2.0 
## [67] codetools_0.2-16    glue_1.3.1          evaluate_0.13      
## [70] lsei_1.2-0          metap_1.1           data.table_1.12.2  
## [73] png_0.1-7           Rdpack_0.11-0       gtable_0.3.0       
## [76] RANN_2.6.1          purrr_0.3.2         tidyr_0.8.3        
## [79] future_1.12.0       assertthat_0.2.1    ggplot2_3.1.1      
## [82] xfun_0.6            rsvd_1.0.0          survival_2.43-3    
## [85] viridisLite_0.3.0   tibble_2.1.1        cluster_2.0.8      
## [88] globals_0.12.4      fitdistrplus_1.0-14 ROCR_1.0-7